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PREFACE 


This project is under the supervision of Professor Ivan I. Mueller, 
Department of Geodetic Science and Surveying, The Ohio State University. 
The Science Advisor is Dr. David E. Smith, Code 921, Geodynamics Branch, and 
the Technical Officer through June 30, 1986, was Mr. Jean Welker, Code 903, 
Technology Applications Center. From July 1, 1986, the Technical Officer is 
Dr. Gilbert D. Mead, Code 601, Crustal Dynamics Project, Space and Earth 
Sciences Directorate. All addresses are at Goddard Space Flight Center, 
Greenbelt, Maryland 20771. 

Because of the change in the contract administration, from now on the 
semiannual reports will cover the periods January-June and July-December, 
instead of October-March and April-September. 
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1. CURRENT TECHNICAL OBJECTIVES 


1. Intercomparison of Different Space Geodetic Measurement and Computational 
Data 

2. Regional Deformations and Relative Plate Motions from Robust Estimation 
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2. ACTIVITIES 


2.1 Earth Rotation Parameter Determination from Different 
Space Geodetic Systems 

This work has been nearly completed during this reporting period. After 
some further development of utility program software for analyzing final 
results, the main simulation experiments were performed, and results and 
conclusions compiled. Although awaiting final approval, the PhD dissertation 
[Archinal, 1987a] has been completed and describes the majority of the work 
and results of this study. Some additional work will probably be done, 

including performing simulation experiments using different orbital and 
observational weighting, and more normal observational data rates. It is 

planned that the just mentioned dissertation will be revised to reflect this 
additional work and published as [Archinal, 1987b], 

During the latter part of the reporting period the results of the study 
were also presented at several meetings, as listed at the end of this report. 


References 


Archinal, Brent A. (1987a), Determination of Earth Rotation by the Combination 
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Univ., Columbus, in preparation. 

Archinal, Brent A. (1987b), "Determination of Earth Rotation by the Combination 
of Data from Different Space Geodetic Systems," Dept, of Geodetic Science and 
Surveying Rep., Ohio State Univ., Columbus, in preparation. 
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2.2 Utilization of Range-Difference Observations 
in Geodynamics 

As anticipated in the 16th Semiannual Report, this study is close to being 
completed, and the final report is in preparation [Dedes, 1987]. For this 
reason, only a brief summary of the work accomplished in the last six months 
will be reported. 

During this period the dynamic models employed by the GEOSPP software 
have been updated to be compatible with the MERIT Standards [Melbourne et 
al., 1983]. GEOSPP is a software used to implement the Simultaneous 
Range-Difference (SRD) semidynamic mode method [Pavlis, 1982], This method 
was anticipated to yield baseline estimates whose accuracy would be 
insensitive to the overall orbital accuracy and sensitive to the accuracy of the 
observations, especially if certain geometric configurations are fulfilled. This 
in fact has proven to be the case. For instance, baselines close to being 
parallel with their observed arcs have been estimated with an accuracy of 2 
cm. This accuracy not only represents one standard deviation but also 

reflects the agreement when these baselines are compared with the 
corresponding ones as estimated by other computational centers using the 
range dynamic mode method [Tapley et al., 1985]. This accuracy was achieved, 
despite the fact that the input initial state vectors were in error by 150 
meters and that the standard deviations of the adjusted ones were as large as 
50 meters. 

These results demonstrate the sensitivity of the SRD method to the 
accuracy of the observations as well as the insensitivity of this method to the 
overall orbital accuracy. Furthermore, these baseline estimates were obtained 
by using only eight Lageos passes. Realistically, eight Lageos passes may be 
observed in about a week. Therefore, this method is ideal for rapid baseline 
determination especially when mobile laser systems are moved into an area 
with the a im of sp endi ng t he le ast possible time collecting data for baseline 
determination. 

There still needs to be established how the relative orientation between 
the baselines and their observed arcs is related to the number of arcs 
necessary for convergence, when the baselines are estimated through the SRD 
method. This will be revealed by estimating baselines having different lengths 
and different orientation with their observed arcs. Furthermore, in the 
geometric mode method (16th Semiannual Status Report) tidal effects should be 
accounted for. Both of these tasks are expected to be completed by the end 
of this calendar year. 

Presently the work is focused on finishing with these tasks and on 
preparing the final report which is expected to be ready within the next three 
months. 
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2.3 An Algorithm for Crustal Deformation Analysis 


Sequential Model Discrimination 


In the previous five Semiannual Reports (12 through 16), a general 
algorithm for the analysis of crustal movement measurements was proposed. 
In this report, a method based on the Bayesian philosophy and entropy 
measure of information is given for the elucidation of time-dependent models 
of crustal motions as part of the proposed algorithm. The method, which is 
due to Box and Hill (1967), is first reproduced; then the strategy of model 
discrimination and design of measurements is illustrated in an example for the 
case of crustal deformation models. 

In the analysis of statistical decision theory two problems are generally 
distinguished. One is the problem of making the best decision on the basis of 
a given set of data, and the other is the problem of designing the best 
experiments in order to get information upon which a decision will be made. 
In this respect, model discrimination procedures of the proposed algorithm can 
be examined within the scope of these two problems. 

Assume that geodetic surveys are performed at different epochs. The 
data indicate that crustal motions occurred during these time intervals and a 
set of concurring descriptive models are postulated either as a result of the 
displayed dislocation patterns of the network points and/or previous 
information. The question is how to use the data to select the best model and 
to design new optimal observations to facilitate model discrimination. 

Since in general every decision is a result of some type of decision rule, 
it is necessary at this point to define a discrimination criterion. 

Consider a complete set of events, E if i = 1, 2, ... m whose probabilities 

m 

are p lf * . . . P m such that £ p i = 1. The expected information of the 

i = l 

message on the occurrence of one of these events is defined as 

I: = 2 Pi In — 7 - ( 1 ) 

1 = 1 Pi 

which is also known as the entropy of the distribution whose probabilities are 
p a i=l, 2, ... m (Shannon 1948). The least possible information occurs when 
Pi = P 2 ••• P m = 1/m, which can be derived by maximizing the above equation 
subject to the constraint £ p 1 = 1# In this case the amount of information is 
small and the entropy as a measure of disorder is maximum. In other words, 
all events are equally likely. In situations where the probability of one event 
Pi is larger than the probability of other events p<, j*i, the amount of 
information is considered to be large and entropy is small. 

This concept can be applied to the discrimination of different competing 
models. Let there be a set of m competing models and the a priori probability 
of the i th model being true is If the observations are performed and the 

a posteriori probability for the i th model is computed, then the information 
gained by this experiment is specified, from ( 1 ) as 
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In Pin-J 


(2) 


AI(H,x) = - X 

i = i 


+ l Pm In p 1n 
1 = 1 


The maximum of AI(H,x) is of interest in order to obtain the greatest amount of 
information out of the experiment. Now, using (2) and considering that there 
exists a finite number of models, the expected change in entropy (information) 
E(AI) = :AJ before and after the n th observation is 


AJ = X Pt„-i 

1 = 1 


f m 

I 

l k = l 


Pkn In p kn 


“ X Pkn - l In p k 

k = 1 


n— i } • Pi(y'x) dy n 


(3) 


where 


Pin-i Pi(y'x) 

m 

X P< n — l Pl(y'x) 

4=1 


(4) 


Substituting (4) into (3) 


AJ = X Pin— t 
1 = 1 


Pi (yix) 

Pi(ylx) In 

.X Pin— l Pl(y'x) 
1=1 


dy« 


(5) 


Now an observation which maximizes (5) is considered the optimum one. 
Evaluation of (5) however is quite complicated, but an upper bound for this 
expression leads to a tractable form. Consider 


m 

X Pjn— l 
1=1 J 


p 4 (ylx) In 


Pi (y<x ) 
Pj(y*x) 


> Pi(y'x) 


In 


Pi (y’x) 

m 

X Pin Pl(y'x) 

< = i 


( 6 ) 


which is a result of Corollary 3.1 of Kullback (1959). Substitution of (6) into 
(5) gives an upper bound AJ U for AJ 


m m 

AJ U = X XPin-l P in— i 

i=i j=l J 


, , Pi(y'x) 

p,<ylx)1 ” JJcTno dy " + 


p j (yl,t>ln dy ™ 


(7) 


Let now a group of competing models is given by 


E(y(0) = a( , )x( < ) , i = 1, 2, ... m (7a) 

where y^ 1 ) is the mxl vector of observations for the i th model, A^') is the 
nonstochastic n*u design matrix, and x^ 1 ) is the uxl unknown parameter 
vector for the i th model, u is not necessarily the same for all models. If the 
observations y n are assumed to be distributed normally with mean E(y„) and 
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known variance <r a , the following relationships hold 


Pi (yJECYn).*) = y== ex p{~ 2^2 [yn - E(y n )] } (8) 

Pi (E(y„) »<r) = j== exp{- y-1 - [E(y n ) - y^ 1 ^] } (9) 

where y^(0 is the predicted value of y n under model i using n-1 observations 
and its variance is <r ( J which is given by 

<r,? = <7 2 a<j0(A(0'A(0) _1 a £0' (10) 

where 1 ) is the row vector of a( f ). From the definition of the probability 
density function of y n under model i given <r and n-1 observations 

Pi(yn'<0 = J Pi (yJECyn) >v) p<(E(y„)l<r) d(E(y n )) (11) 


If (8) and (9) are substituted in (11) and integrated, then 

P < (y " ,<0 = 2n(^ ' V<r " f) exp {" (y " * (12) 

Substitution of (12) into (7) results in the following operational form of 
discrimination function 


\ m m 

AJu - p Z E Pin—i Pjn~i 

* 1=i j=i J 


(<r? - <r1)‘ 


(<r a 


+ <r* ) (<r a 


•1) 


W> - yfJ>) 


f 1 

x } 

U<r 2 + v?) 

(* a + 


(13) 


where y^j 1 ) and y^l) are the predicted observations for model i and j, <j s 2 and 
<Tj a are the predicted variances obtained from (10) for these predicted 
observations and <r 2 is the known a priori variance of the observations. 

It is now possible to design sequential geodetic surveys to discriminate 
the descriptive models of deformations using this entropy measure of 
information. The scenario which is depicted in Fig. 1 is as follows. 

First an initial network design for the area under consideration is 
constructed, for instance, using the D-optimal design criteria which is 
discussed in the 12th Semiannual Report. Geodetic surveys are then 

performed at two different epochs covering the whole network. This is 

followed by the estimation of deformation parameters from the differences of 
observed quantities. At this point, information provided by the current 
estimates and prior qualitative information are examined, and prior 
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probabilities are assigned to each model. If no preference is inferred from 
the existing information, each model is assigned equal probabilities p = 1/m, 
where m is the number of models. The next optimal observation (not 
necessarily the resurvey of the whole network) that gives the maximum 
expected discrimination among m rival models is sought in (13). Then the new 
optimal measurement (s) is performed and posterior probabilities for each model 
are computed using equation (4). Finally, the current standing of each model 
is examined. This procedure is repeated each time using posterior 
probabilities of previous observations as prior probabilities for the succeeding 
observations, until one model emerges from the others. The following 
numerical example illustrates the procedure. 

Numerical Example 

In the previous section, a method based on the entropy measure of 
information is presented as a possible candidate for the discrimination of 
several competing time-dependent models. In order to get a better feeling for 
the applicability of the method to the crustal deformation analysis, this section 
describes a numerical example. 

Consider the case of a homogeneous deformation field. The following 
descriptive models are postulated as a result of prior experiments to represent 
possible network deformations, 


model 

i. 

dx = 
dy = 

(e x x + ©xyY) 

( a xyX a yY ) 

At 

At 

(14) 

model 

2. 

dx = 

e x x At 


(15) 

model 

3. 

dx = 
dy = 

e xa Y At 
e X yX At 


(16) 


where dx and dy are the displacement components of network points for the 
period At; e x and e y are the extensional strains in X Y directions; e x „ is the 
shearing strain. If the baselines are observed at different epochs, Cnen the 
baseline length t j j at epoch t is given by the following expression, 

*ij = (*j - *i) 2 + (Yj ~ Yi) 2 (17) 

Linearizing this expression about the initial epoch t„ and considering 
equations (14), (15) and (16), results in the following mathematical models, 


model 

1: 


- *ijt 0 = At t 

ijCsin^j e x +cos 2 « 

model 

2: 


- **jt 0 = At 

sin 2 otj j e x 

model 

3: 

‘<jt 

" *<1*0 = At 

sin 2a j j e X y 


e a +sin 2a i j 


e xiJ ) 


(18) 

(19) 

( 20 ) 


where *jjt and < ijt 0 are the observed baseline lengths at epochs t and t 0 re- 
spectively, and otjj is the azimuth of the observed baseline i-j. In this 
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example, model 3 is chosen to be the correct model and the pseudo- 
observations and 4 ijt 0 are derived using this model. They are also con- 

taminated with a noise from a normal distribution with mean zero and variance 
1 mm. Shearing strain, e x ^, in the correct model 3 is 0.50 ppm and the time 
interval between observations is constant and equal to a month. The 

sequential model discrimination procedure can now be performed using the 
algorithm depicted in Fig. 1. 

An initial design is set up for the measurement of deformation parameters 
(Fig. 2a). This is a D-optimal design of model 1, which was derived in the 
12th Semiannual Report, with some additional observations. It is initially 
postulated that all three models are equally likely. In other words, prior 
probabilities for each model are 1/3. Part of the network is resurveyed a 
month later. This makes it possible to predict a new observation for the next 
month that maximizes the discrimination function AJ U for each model and select 
the possible observation which gives the maximum discrimination. In this 
example three alternatives are possible: baselines which are in the north-east 

direction, baselines which are in the south-west direction, and baselines which 
are in the east-west direction. The predicted observation is then performed 
and posterior probabilities for each model are computed using the new 
observations and prior probabilities (4). 

The rest of the experiment continues following this prediction and 
observation procedure until computed posterior probabilities indicate that one 
model is superior to the others (Figs. 2c, 2d and 2e). Fig. 2 shows that the 
correct model 3 is identified effectively with 14 baseline observations after six 
months. 

As a comparison, null-hypothesis tests (H 0 : Model 3 is the same as model 

1, H 0 : Model 3 is the same as model 2) are performed using the estimated 

parameters obtained by the least squares method at 5% and 1% significance 
levels. Neither hypothesis is rejected at the 1% significance level until At = 5 
months (Fig. 2e). In the case of <x - 0.05, model 2 is rejected at At = 5. 
Model 2 is rejected at At = 6 months at both levels. However, the results 
were ambiguous in the sense that model 2 and model 3 were still likely 
candidates until the last measurement was performed. This problem is clearly 
eliminated by the proposed method due to the history of accumulated 
measurements and calculated posterior probabilities of each model. 
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